Stability of Numerical Solution of Kinematic Wave Model

Theoretical background, visuals, and step-by-step solutions to numerical problems.

1. Theoretical Background

The Kinematic Wave Model simplifies the 1D Saint-Venant equations by assuming that the friction slope is approximately equal to the bed slope ($S_f \approx S_0$). This means that gravity and friction forces balance each other, and pressure/inertial terms are neglected.

The governing continuity equation is:

$$\frac{\partial A}{\partial t} + \frac{\partial Q}{\partial x} = q$$

Where $Q$ is a function of flow area $A$ based on Manning's equation. For a channel where hydraulic radius $R \approx y$ (flow depth), the discharge is given by:

$$Q = \frac{1}{n} A R^{2/3} S_0^{1/2} = \left( \frac{\sqrt{S_0}}{n B^{2/3}} \right) A^{5/3} = \alpha A^\beta$$

The speed at which the kinematic wave travels is called the Kinematic Wave Celerity ($c_k$):

$$c_k = \frac{dQ}{dA} = \alpha \beta A^{\beta-1} = \frac{5}{3} \frac{Q}{A} = \frac{5}{3} V$$

Courant Condition Visualized on x-t Grid

x t (i, n) (i+1, n) (i+1, n+1) (i, n+1) Characteristic dx/dt = c_k Δx Δt

For stability, the physical wave must not cross the spatial boundary $\Delta x$ within the time step $\Delta t$.

The Stability Criterion (Courant Condition)

To maintain numerical stability in explicit finite-difference schemes, the Courant-Friedrichs-Lewy (CFL) condition must be satisfied. This ensures that the numerical domain of dependence encompasses the physical domain of dependence. The Courant number ($C_r$) must be less than or equal to 1.

$$C_r = c_k \frac{\Delta t}{\Delta x} \le 1 \quad \implies \quad \Delta t_{max} = \frac{\Delta x}{c_k}$$

2. Solution to Question 1

"A channel having a width of 30 m carries a discharge of 110 m³/s through a section. The bed slope and the Manning's 'n' value are 3% and 0.035 respectively and the hydraulic radius is equal to flow depth. Recommend the maximum time step required for stable solution of kinematic wave routing in this condition, if the value of Δx is considered as 1400 meters."

Given Data:

  • Width ($B$) = $30 \text{ m}$
  • Discharge ($Q$) = $110 \text{ m}^3\text{/s}$
  • Bed slope ($S_0$) = $3\% = 0.03$
  • Manning's coeff ($n$) = $0.035$
  • Hydraulic radius ($R$) = $y$ (flow depth)
  • Grid space ($\Delta x$) = $1400 \text{ m}$

Step-by-Step Calculation:

Step 1: Compute Flow Depth ($y$)

Using Manning's Equation: $Q = \frac{1}{n} A R^{2/3} S_0^{1/2}$

Substitute $A = B \cdot y$ and $R = y$:

$$110 = \frac{1}{0.035} \cdot (30y) \cdot (y)^{2/3} \cdot (0.03)^{1/2}$$ $$110 = \frac{\sqrt{0.03}}{0.035} \cdot 30 \cdot y^{5/3}$$ $$110 = \left( \frac{0.173205}{0.035} \right) \cdot 30 \cdot y^{5/3}$$ $$110 = 148.4615 \cdot y^{5/3}$$ $$y^{5/3} = \frac{110}{148.4615} = 0.74093$$ $$y = (0.74093)^{3/5} = 0.8353 \text{ m}$$

Step 2: Compute Flow Area ($A$) and Velocity ($V$)

$$A = B \cdot y = 30 \cdot 0.8353 = 25.06 \text{ m}^2$$ $$V = \frac{Q}{A} = \frac{110}{25.06} = 4.389 \text{ m/s}$$

Step 3: Compute Kinematic Wave Celerity ($c_k$)

For $R=y$ (equivalent to a very wide channel assumption), $\beta = 5/3$.

$$c_k = \frac{5}{3} V = \frac{5}{3} \cdot 4.389 = 7.316 \text{ m/s}$$

Step 4: Compute Maximum Time Step ($\Delta t_{max}$)

Applying the Courant condition:

$$\Delta t_{max} = \frac{\Delta x}{c_k} = \frac{1400}{7.316} = 191.36 \text{ seconds}$$

3. Solution to Question 2

"A channel with a width of 40 m, bed slope 2% and Manning’s n = 0.03 carries a discharge of 100 m³/s through a section. If Δx is taken as 1500 meters, recommend the maximum time step for stable solution of kinematic wave routing in this condition. Assume hydraulic radius equal to flow depth."

Given Data:

  • Width ($B$) = $40 \text{ m}$
  • Discharge ($Q$) = $100 \text{ m}^3\text{/s}$
  • Bed slope ($S_0$) = $2\% = 0.02$
  • Manning's coeff ($n$) = $0.03$
  • Hydraulic radius ($R$) = $y$ (flow depth)
  • Grid space ($\Delta x$) = $1500 \text{ m}$

Step-by-Step Calculation:

Step 1: Compute Flow Depth ($y$)

Using Manning's Equation: $Q = \frac{1}{n} (B \cdot y) y^{2/3} S_0^{1/2}$

$$100 = \frac{1}{0.03} \cdot (40y) \cdot (y)^{2/3} \cdot (0.02)^{1/2}$$ $$100 = \frac{\sqrt{0.02}}{0.03} \cdot 40 \cdot y^{5/3}$$ $$100 = \left( \frac{0.141421}{0.03} \right) \cdot 40 \cdot y^{5/3}$$ $$100 = 188.5618 \cdot y^{5/3}$$ $$y^{5/3} = \frac{100}{188.5618} = 0.53033$$ $$y = (0.53033)^{3/5} = 0.6844 \text{ m}$$

Step 2: Compute Flow Area ($A$) and Velocity ($V$)

$$A = B \cdot y = 40 \cdot 0.6844 = 27.376 \text{ m}^2$$ $$V = \frac{Q}{A} = \frac{100}{27.376} = 3.653 \text{ m/s}$$

Step 3: Compute Kinematic Wave Celerity ($c_k$)

Using the relationship $c_k = \frac{5}{3} V$:

$$c_k = \frac{5}{3} \cdot 3.653 = 6.088 \text{ m/s}$$

Step 4: Compute Maximum Time Step ($\Delta t_{max}$)

Applying the Courant condition:

$$\Delta t_{max} = \frac{\Delta x}{c_k} = \frac{1500}{6.088} = 246.39 \text{ seconds}$$